from math import fabs, log
from random import uniform

# Plan dzialania:
#
# 1. Wyliczyc wszystkie mozliwe rekurencje.
# 2. Wylosowac parametry.
# 3. Wyliczyc wartosc lambdy. 


class Recurrencer:
    @classmethod
    def get_root(cls, d1, d2):
        if (d1 == d2):
            return 2**(1 / d1)
        
        t = min(d1, d2)
        d1 = max(d1, d2)
        d2 = t

        def f(x):
            return x**(d1) - x**(d1 - d2) - 1

        # f' = d1 * x^(d1 - 1) - (d1 - d2) * x^(d1 - d2 - 1)
        # f' = x^(d1 - d2 - 1) * ( d1 * x^(d2) - (d1 - d2) )
        # f' == 0 <==> (x == 0) or ( ( d1 * x^(d2) - (d1 - d2) ) == 0 )
        # f' == 0 <==> (x == 0) or ( ( x^(d2) == (d1 - d2) / d1 )
        # f' == 0 <==> (x == 0) or ( d2 * ln(x) == ln(1 - (d2 / d1) ) )
        # f' == 0 <==> (x == 0) or ( x == e^(ln(1 - (d2 / d1) ) / d2 ) )
        
        x0 = 2**(log(1 - (d2 / d1), 2) / d2)
        x1 = max(x0, 1.0)

        # print 'd1 = %f, d2 = %f, x0 = %f, x1 = %f' % (d1, d2, x0, x1)

        while (f(x1) < 0):
            x0 = x1
            x1 *= 2

        _max_error = 10**-6
        _max_root = 1.5

        while (x1 - x0 > _max_error):
            assert (x1 - x0 > 0)

            x = (x1 + x0) / 2
            if (f(x) > 0):
                x1 = x
            else:
                x0 = x

            if x0 > _max_root:
                return _max_root

        root = (x1 + x0) / 2

        # print 'fin'

        return root

    @classmethod
    def add_elements(cls, set_size, num_elements):
        def add_element(buckets):
            if len(buckets) == 1:
                return [[buckets[0] + 1]]

            l = [buckets[0] + 1] + list(buckets[1:])

            r = [l]

            for l in add_element(buckets[1:]):
                t = [buckets[0]] + l
                r.append(t)

            return r

        l = []
        for _ in range(num_elements):
            l.append(0)

        all_buckets = set([tuple(l)])

        for _ in range(set_size):
            l = []
            for b in all_buckets:
                l.extend(map(tuple, add_element(b)))
            all_buckets = set(l)

        return all_buckets

    def __compute_in_deltas(self):
        self.__in_deltas = \
            (0,

             self.__alphas[0] + self.__betas[0],

             min(2 * self.__alphas[0], 
                 self.__alphas[1]) + 
             self.__betas[0],

             min(3 * self.__alphas[0], 
                 self.__alphas[0] + self.__alphas[1]) + 
             self.__betas[0],

             min(3 * self.__alphas[0], 
                 self.__alphas[0] + self.__alphas[1], 
                 self.__alphas[2]) + 
             self.__betas[0])

    def __compute_dlt_betas(self):
        deltalized = [self.__betas[0]]
        for i in range(len(self.__betas) - 1):
            deltalized.append(self.__betas[i + 1] - self.__betas[i])
        deltalized.append(1.0 - self.__betas[-1])

        self.__dlt_betas = tuple(deltalized)

    # params should cover elements from 2nd, for xmpl. alpha(2), beta(2)
    # even if those are equal 0. elements not covered are assumed to be 
    # equal 1.
    def __init__(self, params, min_set_size = 3):
        assert(len(params['alpha']) >= 3 and len(params['beta']) >= 1)
        assert(min_set_size >= 3)

        self.__alphas = tuple(params['alpha'])
        self.__betas = tuple(params['beta'])

        self.__min_set_size = min_set_size

        self.__num_alphas = len(self.__alphas)
        self.__num_betas = len(self.__betas)

        self.__compute_dlt_betas()
        self.__compute_in_deltas()

        # it's more efficient to have all frequencies calculated separately
        self.__frequencies = {}

        # delta_alpha == 0 and delta_beta == 0 for |S| >= max_set_size
        self.__max_set_size = max(max(self.__num_alphas, self.__num_betas) + 3,
                                  self.__min_set_size)

        # first bucket is for elements with frequency 1, which we acturally
        # don't count. last bucket is for elements with high frequencies that
        # betas don't cover, for which we assume beta is 1
        sieve = self.add_elements(self.__max_set_size, self.__num_betas + 1)

        for i in reversed(range(self.__min_set_size, self.__max_set_size + 1)):
            sieve = [f for f in sieve if sum(f[:-1]) <= i]
            tmp = []
            for l in sieve:
                l = list(l)
                l[-1] = i - sum(l[:-1])
                tmp.append(tuple(l))
            assert(len(tmp) > i)
            self.__frequencies[i] = tuple(tmp)
        

    def reset(self, params):
        assert(len(params['alpha']) == self.__num_alphas and 
               len(params['beta']) == self.__num_betas)

        self.__alphas = tuple(params['alpha'])
        self.__betas = tuple(params['beta'])
        self.__compute_dlt_betas()
        self.__compute_in_deltas()

    def __alpha_s(self, set_size):
        self.__alpha = 0
        if set_size - 1 > self.__num_alphas:
            self.__alpha = 1
        elif set_size > 1:
            self.__alpha = self.__alphas[set_size - 2]
        return self.__alpha

    def __dlt_alpha_s(self, set_size):
        return self.__alpha_s(set_size) - self.__alpha_s(set_size - 1)

    def __rec_out(self, set_size):
        ret = []

        for frequency in self.__frequencies[set_size]:
            assert(len(frequency) == len(self.__dlt_betas))
            # for frequencies > self.__dlt_betas, delta beta == 0
            fxdb = [f * b for (f, b) in zip(frequency, self.__dlt_betas)]

            if frequency[0] < 3:
                delta = self.__in_deltas[frequency[0]]
            elif frequency[0] == 3 and set_size == 3:
                delta = self.__in_deltas[3]
            elif frequency[0] >= 3 and set_size >= 4:
                delta = self.__in_deltas[4]
            else:
                print 'freq[0]: %d, set_size: %d' % (frequency[0], set_size)
                assert(False)

            ret.append(self.__alpha_s(set_size) + sum(fxdb) + delta)
            
        return ret

    def __rec_in(self, set_size):
        ret = []

        for frequency in self.__frequencies[set_size]:
            assert(len(frequency) == self.__num_betas + 1)
            # for frequencies > self.__dlt_betas, delta beta == 0
            fxb = [f * b for (f, b) in zip(frequency, self.__betas)]

            w = [(i + 1) * frequency[i] for i in range(len(frequency))]

            ret.append(self.__alpha_s(set_size) +
                       sum(fxb) +
                       frequency[-1] +
                       self.__dlt_alpha_s(set_size) * sum(w))


        return ret

    def __rec(self, set_size):
        r1 = self.__rec_out(set_size)
        r2 = self.__rec_in(set_size)

        assert(len(r1) == len(r2))

        return zip(r1, r2)

    def get_max_root(self):
        ret = 1.0

        # print 'a, min: %d, max: %d' % (self.__min_set_size, self.__max_set_size)
        for set_size in range(self.__min_set_size, self.__max_set_size + 1):
            r = [self.get_root(d1, d2) for (d1, d2) in self.__rec(set_size)]
            assert(len(r) > 0)
            # r = [(self.get_root(d1, d2), (d1, d2), f) for ((d1, d2), f) in zip(self.__rec(set_size), 
            #                                                                    self.__frequencies[set_size])]
            # print max(r)[1]
            # print len(self.__frequencies[set_size])
            r.append(ret)
            ret = max(r)
            # print 'set_size: %d, max: %f' % (set_size, log(ret, 2))
        # print 'b'
        # print '----'
        # print ret[1]
        # print '----'

        # print log(ret, 2)

        return log(ret, 2)

# ---------- D E B U G -----------
assert_r = Recurrencer({'alpha': (0.377, 0.754, 0.909, 0.976), 
                        'beta':  (0.399, 0.767, 0.929, 0.985)}, 3)

assert(int(assert_r.get_max_root() * 1000) == 305)
# --------------------------------

# r = Recurrencer({'alpha': (0.004, 0.619, 0.842, 0.951), 
#                  'beta':  (0.553, 0.831, 0.945, 0.993)})

# print r.get_max_root()

# r = Recurrencer({'alpha': (0.002, 0.640, 0.849, 0.950, 1.000), 
#                  'beta':  (0.463, 0.768, 0.925, 0.999, 1.000)})


# print r.get_max_root()

# r = Recurrencer({'alpha': (0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 1.000), 
#                  'beta' : (0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 0.001, 1.000)})


# print r.get_max_root()



# ----------------------------------------------------

def seek_best_params(min_set_size, max_set_size):
    def rand_params(seed, params_range):
        ret = []
        last = 0
        for i in seed:
            ret.append(uniform(max(last, i - params_range),
                               min(i + params_range, 1)))
            last = ret[-1]
        return tuple(ret)

    iterations = 100
    min_range = 0.00001

    range_var = 1.0
    ones = tuple([1 for i in range(max_set_size - 3)])
    ret = {'alpha': ones, 'beta': ones}
    r = Recurrencer(ret, min_set_size)
    ret_root = r.get_max_root()

    while (range_var > min_range):
        print 'yield'
        for i in range(iterations):
            tmp = {'alpha': rand_params(ret['alpha'], range_var), 
                   'beta': rand_params(ret['beta'], range_var)}
            r.reset(tmp)
            tmp_root = r.get_max_root()
            if tmp_root < ret_root:
                ret = tmp
                ret_root = tmp_root
                print ret_root
                print ret
                print range_var
        range_var /= 2

    print ret_root
    return ret

print seek_best_params(8, 8)

    
# from cProfile import run
# run('print seek_best_params(7)')
